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The amount of ionic current flowing through channels is determined by the interplay between two separate 
time-dependent processes: activation and inactivation gating. Activation is concerned with the stimulus-dependent 
opening of the main intracellular gate, whereas inactivation is a spontaneous conformational transition of the 
selectivity filter toward a nonconductive state occurring on a variety of timescales. A recent analysis of multiple 
x-ray structures of open and partially open KcsA channels revealed the mechanism by which movements of the 
inner activation gate, formed by the inner helices from the four subunits of the pore domain, bias the conforma- 
tional changes at the selectivity filter toward a nonconductive inactivated state. This analysis highlighted the important 
role of PhelOS, a residue located along the inner helix, near the hinge position associated with the opening of the 
intracellular gate. In the present study, we use free energy perturbation molecular dynamics simulations (FEP/MD) 
to quantitatively elucidate the thermodynamic basis for the coupling between the intracellular gate and the 
selectivity filter. The results of the FEP/MD calculations are in good agreement with experiments, and further 
analysis of the repulsive, van derWaals dispersive, and electrostatic free energy contributions reveals that the ener- 
getic basis underlying the absence of inactivation in the F103A mutation in KcsA is the absence of the unfavorable 
steric interaction occurring with the large flelOO side chain in a neighboring subunit when the intracellular gate is 
open and the selectivity filter is in a conductive conformation. Macroscopic current analysis shows that the IIOOA 
mutant indeed relieves inactivation in KcsA, but to a lesser extent than the F103A mutant. 



INTRODUCTION 

channels control the flow of ions across cell mem- 
branes principally via two different mechanisms involving 
structural rearrangements either at the intracellular 
gate or at the selectivity filter. The conformational states 
of these two elements are dynamically coupled. Opening 
of the intracellular gate via external stimuli (activation) 
results in a transient period of sustained ion conduc- 
tion until the selectivity filter undergoes a spontaneous 
conformational change (inactivation) toward a collapsed 
or constricted nonconductive state (Hoshi et al., 1991; 
Yellen, 1998; Gao et al., 2005; Cordero-Morales et al., 
2006a, 2007; Cuello et al., 2010b,a). Activation and 
inactivation gating are allosterically coupled: an open 
gate favors an inactivated filter, and a closed gate stabi- 
lizes a conductive filter (Panyi and Deutsch, 2006; Cuello 
et al., 2010a). Ultimately, the rates of entry into and re- 
covery from inactivation govern the total amount of 

membrane permeation available in response to a 
given external stimulus. The complex interplay between 
activation and inactivation, a hallmark of all channels. 
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underlies a wide range of electrical events in excitable 
and nonexcitable cells. Although it is customary to discuss 
activation and inactivation of channels in terms of 
two functional "gates," the molecular basis of how those 
structural elements might be allosterically coupled is not 
well understood. 

The bacterial KcsA K^ channel from Streptomyces lividans 
provides an excellent prototypical model system to 
understand activation-inactivation coupling. Although 
the KcsA channel comprises only a pore domain without 
the complex machinery of voltage-gated K^ channels, it 
contains all the molecular elements necessary for ion 
conduction, activation, and inactivation gating. In the 
closed state, the inner helices of each of the four sub- 
units constituting the pore domain adopt a straight 
configuration, forming a tight bundle that occludes the 
permeation pathway at the intracellular end (Doyle et al., 
1998; Zhou et al., 2001). Opening of the intracellular 
gate is associated with a hinge-bending motion that 
allows these four helices to splay outward, creating a 
wide aqueous vestibular entryway leading up to the narrow 
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selectivity filter (Perozo et al, 1999; Liu et al., 2001; 
Jiang et al., 2002; Cuello et al., 2010b). 

A recent analysis based on multiple x-ray structures of 
open and partially open KcsA revealed the mechanism 
by which movements of the inner bundle activation 
gate is linked to conformational changes within the 
selectivity filter (Cuello et al., 2010a,b). The analysis 
highlighted the important role of PhelOS, a residue 
located along the inner helix near the hinging position 
associated with the opening of the intracellular gate. In 
the open state, it was observed that the aromatic side 
chain of PhelOS changes its rotameric state, leading to 
an apparent steric clash with residues of the selectivity 
filter from the neighboring subunit (Fig. 1). Mutations 
at position 103 were shown to affect gating kinetics in a 
size-dependent way, with small side-chain substitutions 
leading to considerably slower entry into inactivation. 
Strikingly, an alanine residue substituted at position 103 
abolished inactivation in KcsA (Cuello et al., 2010a). 
Average energies calculated from equilibrium molecu- 
lar dynamics (MD) simulations indicated that, among 
the different residues along the inner helix, the strongest 
energetic interactions with the selectivity filter involved 
the side chain of Phel03 (Cuello et al., 2010a). These 
observations led to the proposition that Phel03 is the 
structural link responsible for the allosteric communication 




Figure 1 . (A) Structures of one subunit of KcsA with its intra- 
cellular gate closed (blue) and open (yellow) . The side chain of 
PhelOS is highlighted in a licorice representation. Zoomed-in 
representations of two adjacent subunits of KcsA in the closed (B) 
and open (C) states showing the residues important for allosteric 
communication between the intracellular gate and the selectivity 
filter as spheres. 



between the intracellular gate formed by the inner heli- 
cal bundle and the selectivity filter. 

Although this analysis based on static x-ray structures 
and equilibrium MD trajectories provided a plausible 
and coherent picture, the functional and crystallographic 
results beg for a more quantitative characterization of 
the nature of the microscopic interactions involving the 
residue at position 103 with the rest of the protein struc- 
ture, and how it underlies the allosteric mechanism 
coupling the activation gate to the selectivity filter. In 
particular, knowledge of average interaction energies as 
calculated by MD is insufficient to explain the factors 
underlying the structural stability of functional states. 
The latter is governed by relative free energies, and this 
is the critical information that is ultimately needed for a 
complete and quantitative analysis of the wild-type (WT) 
and mutant channels. 

Here, alchemical free energy perturbation (FEP) MD 
simulations (FEP/MD) of the F103A mutation in KcsA are 
used to elucidate the thermodynamic basis for the cou- 
pling between the intracellular gate and the selectivity 
filter. The results of the calculations are in good agree- 
ment with experiments. This then enables us to carry 
out a detailed dissection of the microscopic interactions 
underlying the relative stability of the two functional 
states of the filter. Decomposition of the free energy 
differences in terms of repulsive, van der Waals, and 
electrostatic contributions reveals the energetic basis 
underlying the absence of inactivation in the F103A 
mutation in KcsA. When the intracellular gate is open, 
the large aromatic phenylalanine side chain at position 
103 appears to destabilize the conductive conformation 
of the selectivity filter because of an unfavorable steric 
van der Waals clash with He 100 in an adjacent subunit. 
This idea is verified experimentally by examining the 
behavior of the lOOA KcsA mutant, thus lending sup- 
port to a mechanical coupling mechanism of activation 
and inactivation that is highly cooperative in nature. 

MATERIALS AND METHODS 

Electrophysiology 

Experimental macroscopic currents through the KcsA channel 
were recorded in response to an activating pH jvimp (from 8 to 3) 
using the liposome patch-clamping methodology that was de- 
scribed previously (Cordero-Morales et al., 2006a). For KcsA-WT 
and F103A mutant channels, macroscopic currents were obtained 
at +100 and —100 mV from KcsA-reconstituted liposomes in 
200 mM of symmetrical conditions. 

Adiabatic energy maps 

Adiabatic energy maps were calculated from scanning the rota- 
meric energy surface along the 5Ci (N-CA-CB-CG) and 
(CA-CB-CG-CDl) dihedral angles of PhelO.3 in one of the four 
subunits of KcsA in its closed (Protein Data Bank accession no. 
1K4C) and open (Protein Data Bank accession no. 3F5W) confor- 
mations. For the sake of simplicity, the interaction was computed 
by including only a single subunit of the KcsA channel in vacuum. 
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Atoms at distances greater than 10 A from the center of mass of 
PhelOS were fixed in space; all other non-hydrogen atoms, except 
for those belonging to Phel03, were harmonically restrained with 
a force constant of 100 kcal/mol/A^. The values of Xi and X2 were 
systematically varied along a grid, and new conformations of the 
PhelOS side chain were built using the internal coordinates module 
of the CHARMM program (Brooks etal., 2009). At eveiy grid point, 
the angles were restrained with a force constant of 1,000 kcal/mol, 
the structure was minimized with 1,000 steps of steepest descent 
and 500 steps of adopted basis Newton-Raphson minimization, 
and the final energy was recorded without the dihedral restraint. 

Atomic model of the KcsA simulation systems 

The x-ray structure of the KcsA K* channel in the open-inactivated 
(oi) state with the intracellular gate open by 32 A, as measured 
from the distance between the Cas of Thrll2 of opposing sub- 
units (Protein Data Bank accession no. 3F5W), was used 
(Cuello et al., 2010b). Missing atoms at the N and C termini of 
3F5W were modeled from the x-ray structure of KcsA in the closed 
state (Protein Data Bank accession no. 1K4C) after local align- 
ment (Zhou etal., 2001). The structure of the channel with an open 
intracellular gate and open-conductive filter (00) was modeled 
by replacing residues 70-80 of the selectivity filter with the co- 
ordinates from the x-ray structure with a closed intracellular gate 
and an open-conductive filter (Protein Data Bank accession no. 
1K4C) (Zhou et al., 2001). In both cases, the simulated model 
comprises the channel structure embedded in a dipalmitoylphos- 
phatidylcholine (DPPC) bilayer surrounded by an aqueous KCl 
salt solution. The models contained the KcsA tetramer (404 
amino acid residues), 112 DPPC molecules, 7,325 water mole- 
cules, and 2 K* ions in the pore at positions SI and S4, and S2 and 
S4, for inactive and conductive filters, respectively. To ensure elec- 
trical neutrality and mimic a 150-mM KCl concentration, 18 K* 
and 32 CP ions were added to the bulk solution. The glutamic 
acid at position 71 was protonated (Berneche and Roux, 2002), 
and the termini were patched with standard N- and C-terminal 
capping groups. For stage 1 of the FEP calculations, an alanine 
residue was substituted at position 103 in all four KcsA subuniLs and 
carefully equilibrated. All simulation systems were set up using 
the CHARMM program (Brooks et al., 2009) according to a meth- 
odology described previously (Berneche and Roux, 2000). 

Alchemical FEP/MD simulations 

The free energy for alchemically mutating an alanine into a phe- 
nylalanine at position 103 was calculated as a step-by-step revers- 
ible work. The alchemical transformation was performed in each 



of the four subunits of the KcsA channel, and thus, the FEP/MD 
computation results in the total observable relative free energy 
caused by a F103A mutation inserted into a tetrameric channel. 
Each step consists of an independent MD simulation correspond- 
ing to a different microscopic state, which allows the computation 
to be run efficiently in parallel. The FEP/MD simulations were 
performed on a fully solvated KcsA system in a lipid bilayer with 
ions as described above. Spurious instabilities can occur when 
simulating a moderate to low-resolution structure, which might 
take a long time to equilibrate, thus causing long-time conver- 
gence problems in the FEP/MD calculations. In view of the 
moderate resolution (3.2 A) of the x-ray structure for the 
open-inactivated state (Cuello et al., 2010b), a weak harmonic 
restraint of 0.5 kcal/mol/A^ was applied to non-hydrogen atoms 
of the backbone at a distance of 10 A away from the four PheI03 
residues during the FEP/MD calculations to avoid deviations 
away from the open-channel structure determined by x-ray 
crystallography (Cuello et al., 2010b). This is justified in the 
present case because we aim to extract the most meaningful 
information from the FEP/MD simulations in a reasonable amount 
of computational time. 

The calculation involved two major stages, depicted in Fig. 2. In 
the first stage, the FEP/MD simulations were performed with only 
an alanine at position 103, without the presence of a virtual aromatic 
group attached to the carbon P atom (CB) and decoupled from 
the surroundings. The methyl group (CH3) of the alanine was 
transformed into a methylene group (CH2). The second stage 
involves materializing the remainder of the phenylalanine side 
chain. This part of the calculation was performed with only a phe- 
nylalanine at position 103, without the presence of a virtual HB3 
atom attached to the CB and decoupled from the surroundings. 
In practice, these non-interacting atoms would not affect the cal- 
culated relative free energies and can be left out for the sake of 
simplicity. Because the same FEP/MD procedure was followed for 
both the 00 (open inner gate with open-conductive filter) and oi 
(open inner gate with inactivated filter) states of the channel 
leading to appropriate endpoints, the procedure leads to a prop- 
erly closed thermodynamic cycle. In step la, a harmonic dihedral 
restraint with a force constant of 20 kcal/mol/ rad^ was applied on 
the Xi of alal03. In step lb, the atom type of CB was changed 
from a methyl group carbon (CT3) to that of a methylene group 
atom (CT2), and the Lennard-Jones (LJ) 6-12 interactions of 
HB3 were switched off as its charge was transferred to the CB in 
the presence of the Xi dihedral restraint. 

In steps 2a, 2b, and 2c, the repulsive, van der Waals dispersive, 
and electrostatic contributions to the nonbonded interactions of 



Step la 




Step 1 b 




Step 2d 



Figure 2. Overall alchemical FEP/MD strategy for transforming an alanine into a phenylalanine. Transparent atoms do not interact 
with other atoms in the system. Red dashed lines indicate a set of four atoms whose dihedral angle is harmonically restrained. Details of 
the calculation are described in the section Alchemical FEP/MD simulation. 
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the aromatic group with its surroundings were progressively 
switched on in the presence of harmonic dihedral restraint po- 
tentials acting on the Xi (N-CA-CB-CG) and (CA-CB-CG-CDl) 
dihedral angles of the phenylalanine side chain. 

The dihedral restraint was applied with a force constant of 
20 kcal/ mol/ rad^, using the x-ray structure as a reference. For steps 
2a and 2b, the LJ 6-12 potential was separated into purely repulsive 
and dispersive components using the Weeks-Chandler-Andersen 
(WCA) prescription (Weeks et al., 1971; see also Deng and Roux, 
2004). Specifically, the core repulsion diverges as r''^ at short 
distance, whereas the van der Waals dispersion interaction goes as 
r'^ at long distance. The WCA separation of the LJ 6-12 potential 
is a convenient method to construct these two components based 
on the first derivative (i.e., force), which simplifies the interpreta- 
tion of the free energy decomposition (Deng and Roux, 2004). 
Lastly, all dihedral restraint potentials were removed in the final 
step, 2d. Such dihedral restraints help improve the sampling and 
convergence by reducing the amount of fluctuations during 
intermediate stages of the FEP/MD calculation for steps 2a, 2b, 
and 2c, and do not affect the final relative free energies. 

The progression of the repulsive, dispersive, and electrostatic 
contributions in the oo and oi states (steps 2a, 2b, and 2c) are shown 
in Fig. 5. For each stage, two independent simulations were run 
for the initial and final states of each window, and the total simula- 
tion time was 100 ps. Data from the first 10 ps was discarded for 
eqviilibration. All the FEP/MD calculations were performed using 
the PERT module of the CHARMM program (Brooks et al., 2009). 
The final results are given in Table II, with a complete breakdown 
of the various contributions to the total free energy. 

RESULTS AND DISCUSSION 

The recent crystallographic data on KcsA shows that the 
opening of the intracellular gate is coupled to a confor- 
mational change taking place within the selectivity filter 
toward a constricted nonconductive state (Cuello et al., 
2010b). Furthermore, a change in the rotameric state 
of the aromatic side chain of PhelOS is also observed. 



The closed (1K4C) and fully open channel (3F5W) 
structures with the highlighted PhelOS are shown in Fig. 1 . 
By virtue of its location, PhelOS is an obvious candidate 
residue that could relay information about the confor- 
mation of the intracellular gate to the selectivity filter. 
It is located along the inner helix, in the vicinity of 
both the selectivity filter and the helix hinge bending 
associated with the opening of the intracellular gate. It 
is unclear, however, why the rotameric state of PhelOS 
changes upon the opening of the activation gate. 

To gain more insight into the microscopic factors re- 
sponsible for the change in the rotameric state of PhelOS, 
we compared the two-dimensional maps of the energy 
surfaces along the Xi. Xi dihedral angles of the PhelOS 
side chain in the closed and open state of KcsA. We em- 
phasize that, although these are only adiabatic (energy- 
minimized) maps and not free energy surfaces, the results 
are expected to be meaningful and informative because 
the allowed conformations of the side chain are domi- 
nated by local steric effects. Fig. S shows that that there 
is a shift in the position of the energy minimum of the 
side chain when the intracellular gate opens. For the 
purpose of this study, the open-inactivated x-ray struc- 
ture of the KcsA channel with its intracellular gate 
opened by S2 A (Protein Data Bank accession no. SF5W) 
is used. When the intracellular gate opens, the energeti- 
cally favorable rotamer of the PhelOS side chain changes 
from an "up" state (Fig. S, A and B, white square: Xi = 
— 68; X2 = 105), in which the aromatic ring is tucked 
along the inner helix in close proximity to the selectiv- 
ity filter, to a "down" state (white triangle: Xi = —155; 
Xa = 41), in which the aromatic ring points away from 
the helix. The change reflects the preferred rotamers 




-180 -130 



Figure 3. Adiabatic energy-minimized Xi' Xa maps of PheIO.S in KcsA with the intracellular gate (A) closed (1K4C) and (B) open 
(Protein Data Bank accession no. 3F5W). The maps are calculated by including only the interactions within a single monomer. The 
white square and triangle indicate the rotameric state of F103 in the crystal structures of the closed (1K4C) and open (3F5W) channels, 
respectively. (C) A snapshot of KcsA with the intracellular gate open corresponding to the white square in B. The residues that clash 
with FI03 are highlighted with spheres. 
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commonly seen in the crystal structures of the closed 
(1K4C) and open (3F5W) states (see Fig. 1 A). An in- 
teresting exception, however, is the synthetic KcsAwith 
a D-alanine substituted at the position of Gly77 (Protein 
Data Bank accession no. 2IH1). This is the single oc- 
currence of a crystal structure corresponding to a KcsA 
with its intracellular gate closed and an open-conductive 
selectivity filter with the ring of PhelOS pointing down 
(Valiyaveetil et al., 2006) . This suggests that both rotam- 
ers of PhelOS are allowed when the intracellular gate 
is closed, although the up state appears to be more fre- 
quently observed. This agrees with the computed adiabatic 
map of KcsA with its intracellular gate closed (Fig. 3 A) ; 
although both the up and down states of F103 are ener- 
getically allowed, the up state is more favorable. 

The energy maps confirm that these two rotameric 
states are favored in their respective conformations and 
that the down state is highly favored when the intracel- 
lular gate is open. The change in rotameric state is es- 
sentially driven by the conformation of the backbone 
associated with channel gating. In the two-dimensional 
maps, the up rotameric state observed in the closed 
channel (Fig. 3 A, white square) corresponds to an un- 
favorable energy region in the open channel (Fig. 3 B, 
white square). As shown in Table I, the distance be- 
tween the bulky residues of the inner helix (He 100 and 
Phel03) and those located near the tip of the selectivity 
filter (e.g., Ala73 and Thr74) decreases significantly, by 
^1.5 to 2.0 A, in going from the closed state to the 
open-inactivated state. Paradoxically, although the inner 
helices bend and splay apart upon the opening of the 
intracellular activation gate, there is actually a slight com- 
pression of the upper region of the inner helix (above 
the sharp bend) against the bottom of the selectivity fil- 
ter. This is shown in Table I; there is a reduction of the 
distance between the bulky residues He 100 and Phel03 
of the inner helix and Ala73, Thr74, and Thr75 at the 
bottom of the selectivity filter. These backbone displace- 
ments have a strong impact on the accessible rotamers 
of Phel03. Examination of the conformations corre- 
sponding to the high energy region in Fig. 3 shows that 
the side chain of Phel03 cannot remain in the same 
rotameric state ( Xi = —68; Za = 105) because of clashes 
with Ala73 and Thr74 (junction between the pore helix 
and selectivity filter) , and Leu36 (outer helix) (Fig. 3 C) . 
In the closed state (1K4C), the Ca-Ca distance between 
Phel03 and Ala73 is 9.2 A, and it is 8.3 A between 
Phel03 and Thr74. In contrast, those distances are re- 
duced to 6.8 and 6.9 A in the open state. As a result, 
the rotameric state of Phel03 in the closed state becomes 
energetically forbidden due to the steric clash between 
the aromatic side chain and some residues within the 
same subunit. Thus, in response to outward bending 
of the inner helix associated with the opening of the 
intracellular activation gate, the side chain of Phel03 
is essentially obligated to change its rotameric state. 



By virtue of its response to the intracellular gate and 
its spatial proximity to the selectivity filter, Phel03 
represents an important structural link that couples 
channel opening and C-type inactivation. To transmit 
information from the intracellular gate to the filter, 
however, the change in the Phel03 ro tamer needs to 
affect the stability of the conductive pore. To help 
formulate this issue more quantitatively, it is useful to 
consider the minimal number of channel states that 
could possibly result from these two dynamical elements. 
The scheme involves transitions between four states: 
closed inner gate with open-conductive filter (co) , open 
inner gate with open-conductive filter (oo) , open inner 
gate with inactivated filter (oi), and closed inner gate 
with inactivated filter (ci). In this way, the observed 
kinetic behavior of WT and mutant KcsA can be de- 
scribed in terms of equilibrium population shifts among 
the four preexisting states (co, oo, oi, and ci) (Yifrach 
et al., 2009; Cuello et al., 2010a). The effectiveness of 
the mutation to alter entry into inactivation of the KcsA 
channel once the intracellular gate is open is governed 
by the probabilities P(oo) and P(oi). It is advantageous to 
address the issue from a thermodynamic standpoint by 
considering the relative free energy of the oo and oi states. 



^'^(oo),(oi) ~ ^(oo) ^(oi) 



-k^Tln 



p 



(1) 



Such an expression also holds for mutants of the 
channel. 

The experimentally observed effects of substituting 
the phenylalanine at position 103 by an alanine observed 
experimentally are shown in Fig. 4. The macroscopic 
current response upon a pH jump shows that the WT 
channel inactivates fairly rapidly, whereas the mutant 

TABLE I 
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F103A remains predominantly conductive. The muta- 
tion F103A sharply reduces the rate and extent of C-type 
inactivation. Qualitatively, inactivation appears to be 
essentially abolished by the F103A mutation. However, 
it is necessary to characterize this observation more 
quantitatively to compare with the results of the free 
energy simulations. The most reliable measure of the 
equilibrium effect of the mutation is via the 7/ /max ratio. 
This allows us to estimate P(oo), the steady-state probabil- 
ity of the channel to have its filter in the conductive 
state while the inner activation gate is open. For the 
sake of simplicity, we consider the behavior of the chan- 
nel under a positive applied potential (outward currents) 
in the present analysis. As observed in Fig. 4, the inacti- 
vation of the KcsA selectivity filter is partly voltage 
dependent. At negative potential, the channel inactivates 
more rapidly because there are additional processes 
associated with the side chain Glu71 (Cordero-Morales 
et al., 2006b). It is sufficient at this point to limit the 
analysis to the outward currents. Based on Fig. 4, the se- 
lectivity filter of KcsA-WT has a probability of '^0.05 on 
average to be in the conductive state when the inner 
activation gate is open; P(oo) ranges from 0.02 to 0.07 
based on the analysis of three traces (not depicted) . In 
contrast, the selectivity filter of the F103A KcsA has a 
probability of ^0.90 to be in the conductive state when 
the inner activation gate is open; P(oo) ranges from 0.86 
to 0.90 based on the analysis of five traces (not de- 
picted). For the WT channel, these results imply that 
the probability ratio of the oo state relative to the oi 
state is '^O. 05/0. 95 at positive voltage, corresponding to 
a AG(oo),(oi) that is roughly equal to +1.8 kcal/mol for the 
WT channel. For the F103A mutant, the probability ratio 
of the oo state relative to the oi state increases to 
'^0.90/0.10, yielding a free energy AG(oo),(oi) of roughly 
— 1.3 kcal/mol. Combining the results from several 
macroscopic current traces for WT KcsA and the F103A 
mutant (not depicted) , the net effect of the mutation 
on the relative stability of the oo and oi state, AAG, 
is found to vary between —3.6 and —2.6 kcal/mol 



(the number is negative because the mutation makes 
the state oo more stable relative to the state oi). 

Computations can be used to gain more insight into 
the role of the F103 side chain and the coupling be- 
tween the activation gate and the selectivity filter. A 
simple "brute force" MD simulation approach would 
seek to generate long trajectories to observe directly the 
relative probabilities P(oo) and P(oi) for the WT and mu- 
tant channel. Because inactivation in KcsA is a slow 
process taking place over hundreds of milliseconds, how- 
ever, such an approach would not be practical. Alter- 
natively, the problem can be treated thermodynamically 
by considering the relative free energy of the different 
functional states. Accordingly, the observable effect of a 
mutation on the relative probabilities of the oo and oi 
state can be expressed as. 



(2) 



AG, 



FlOSA.WT 
(oo) 



AG 



F103A,WT 
■(oi) 



Quantities such as ag^^ 



F103A 
'(oo) 



;] and 



^gFi()3A,wT ^ ^c,™^^ _ G™] represent free energy differences 
between the WT and mutant channel in a given confor- 
mational state (oo or oi) . 

Although experimentally inaccessible, such free 
energy differences can be calculated via alchemical 
FEP/MD simulations, in which the step-by-step revers- 
ible thermodynamic work for mutating an amino acid is 
calculated. Each stage consists of an independent MD 
simulation with a different set of options and parameters. 
The physically intuitive route would compare the change 
in the probabilities for inactivation in the WT and mutant 
(F103A). The alchemical FEP/MD simulation approach 
involves comparing the difference in free energy for 
mutating (in silico) the aJanine to a phenylalanine in the 
inactive and open-conductive states. 

The results of the FEP/MD simulations are given in 
Table II. The total free energy difference to replace the 
phenylalanine at position 103 by an alanine in KcsA is 



Figure 4. Macroscopic pH-de- 
pendent currents traces (15 s) 
obtained at 100 mV (A) and 




— 100 mV (B) from proteolipo- 
somes containing KcsA-WT (black) 
and KcsA-F103A mutant (blue). 
WT KcsA is rapidly inactivating, 
whereas the mutant F103A is 
effectively non-inactivating. The 
probability P(oo) for KcsA-WT 
ranges from 0.02 to 0.07 (based 
on the analysis of the ///max ratio 
for three traces), and the prob- 
ability for the F103A KcsA 
mutant ranges from 0.86 to 0.90 
(based on the analysis of the 
///max ratio for five traces) . 
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consistent with a stabilization of the conductive state 
relative to the inactive state by —5.3 kcal/mol. This 
result is consistent with the important observation that 
the F103A mutation relieves inactivation (Cuello et al., 
2010a). Moreover, the calculated AAG is also in good 
accord with the relative free energy estimated from 
macroscopic current experiments, which ranges from 
— 3.6 to —2.6 kcal/mol (see above). Once it is estab- 
lished that the computational result for the F103A 
mutation is consistent with experiment, the FEP/MD 
calculations give us an opportunity to gain a deeper in- 
sight into the molecular interactions driving channel inac- 
tivation. This can be achieved by considering the various 
free energy components of the step-by-step alchemical 
F103A transformation used in the FEP/MD simulations. 

A particularly useful framework for dissecting the 
interaction of a side chain with its surroundings in 
terms of repulsion, dispersion, and electrostatics can be 
formulated based on methods developed for the cal- 
culation of the solvation of molecular solutes (Deng 
and Roux, 2004) . The decomposition proceeds from 
the general observation that all intermolecular forces 
are dominated by harsh short-range repulsive interac- 
tions arising from the Pauli exclusion principle, long- 
range van der Waals attraction as a result of quantum 
mechanical dispersion, and electrostatic interactions 
from a non-uniform molecular charge distribution. In 
standard biomolecular potential functions, electrostatic 
interactions are represented on the basis of Coulomb 
interactions between partial charges, whereas nonpolar 
interactions are modeled with LJ 6-12 potentials. The 
latter can be separated into purely repulsive and attrac- 
tive parts according to the WCA prescription (Weeks 
et al., 1971; Deng and Roux, 2004). In the following, the 
repulsive and dispersive contributions are based on the 
WCA separation of the LJ potential. If the repulsive, 
dispersive, and electrostatic parts of the interactions 
between the side chain and its surroundings are ther- 
modynamically switched on in a step-by-step fashion 
using an FEP/MD simulation, the free energy difference 
arising from the nonbonded interactions may be writ- 
ten as, AAGi„ = AAG„p + AAG,,, + AAG,i„, where 

Grep=LAG(„„) (rep)-AG(„i) (rep)J, 



and 



AAG., 



:[AG, 



F103A,WT/ , X A/^F103A, 

(elec) - AG, 



(oo) 



■^(elec)]. 



Such a free energy decomposition must be inter- 
preted with caution, however, as the various terms are 
conditionally defined by the order in which the interac- 
tions are introduced in the FEP/MD simulations (see 
Materials and methods). Fundamentally, any thermo- 
dynamic decomposition of a free energy in a strongly 
coupled system is a theoretical construct that is path 
dependent. The step-by-step FEP/MD method that we use 
here consists in turning on the repulsive interaction, 
then the van der Waals dispersion, and then the electro- 
statics. This path makes it possible to attribute the free 
energy contributions to a specific interaction, but it is 
dependent on the sequence chosen. For example, the 
dispersion is calculated once the repulsion is already 
active, and the electrostatics is calculated once both the 
repulsion and the dispersion are already active. 
Changing the order of how these interactions are 
turned on, however, could lead to nonsensical results. 
For example, it is not possible to turn on the dispersion 
or the electrostatics before the core repulsion is already 
active because the energy becomes divergent. More 
importantly, the step-by-step FEP/MD decomposition 
method can provide useful insight, as demonstrated by 
previous applications to the solvation free energy of 
simple compounds (Deng and Roux, 2004) , and on the 
binding of ligands to proteins (Deng and Roux, 2009). 

Table II gives the results of the FEP/MD calculations 
with a complete breakdown of the various contributions 
to the total free energy. It is observed that the domi- 
nant contribution arises from the repulsive component 
from the aromatic side chain ( — 3.4 kcal/mol) . As shown 
in Fig. 5 and in Table II, the repulsive component is 
actually larger when the selectivity filter is in its open-con- 
ductive state compared with the inactivated state. As a 
result, the AAGfor the F103A mutant favors the conduc- 
tive state because of the larger repulsive free energy 
interaction of the side chain with the surroundings. 

Thus, the final FEP/MD results are dominated by the 
unfavorable repulsive free energy of the active state rel- 
ative to the inactive state: 



AAG, 



:[AG, 



F103A,WT 
(oo) 



(dis) 



A ^F103A,WT 
■ A^oi) 



(dis)]. 



AAG 



J 0 



aAG™^^5) aAG,^(5) 



(3) 



TABLE II 

Breakdown of the vario us steps in the FEP/MD 



Free energy 


la 


lb 


2a 


2b 


2c 


2(1 


Total 




(Dihedral) 


(Hydrogen) 


(Repulsion) 


(Di.spersion) 


(Electrostatics) 


(Dihedral) 




AG„i 


.S.92 


-3.87 


39.4 


-56.8 


-22.7 


5.94 


-.S8.0 




.S.99 


-3.07 


42.8 


-55.6 


-22.7 


5.95 


-.S2.6 


AAG 


0.07 


-0.76 


-.S.40 


-1.20 


0.00 


0.01 


-5.3 
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Conductive (First 40 ps) 

Conductive (Final 50 ps) 

Inactive (First 40 ps) 

Inactive (Final 50 ps) 



Figure 5. Progression of the energetic contribution (kcal/mol) from the repulsive, dispersive (van der Waals), and electrostatic inter- 
actions based on the FEP/MD simulations as a function of the thermodynamic coupling parameters s, and X. The parameter s serves to 
switch the repulsive contribution, the parameter £, serves to switch the van der Waals dispersion, and the parameter X serves to switch the 
electrostatic interactions. The repulsive and dispersive contributions are based on the WCA separation of the LJ potential. These three 
panels correspond to steps 2a, 2b, and 2c depicted in Fig. 2). The dominant contribution arises from the repulsive interaction. The close 
overlap of the solid and dashed lines, corresponding, respectively, to the first 40 and final 50 ps of the calculation in each window, is an 
indication of the small error in the results. 



where the thermodynamic coupling parameter s serves to 
switch on (i= 1) or off (5=0) the repulsive contribution. 
The largest contribution to the difference in repulsive 
contributions arises when the difference in the derivative, 
[dAG™^'' /ds-dAG^{s)/ds], reaches a maximum. Analysis 
of the configurations from FEP/MD trajectories shows 
that this occurs when the repulsive thermodynamic cou- 
pling parameter, s, is near 0.5. To clarify the origin of 
these contributions, it is useful to rewrite the thermo- 
dynamic derivatives in terms of explicit WCA-repulsive 
LJ interactions, 

AAG^^p = J&j47ir2rfr2^[p<7'(r;5)-p<°"(r;5)], (4) 

0 0 KT 

where ^^yi'"'-^) is the ^-dependent soft-core WCA repul- 
sive potential, and Pai\m) and p',r,!'(r;.s) are the average 
radial distribution functions between atom a (on 
PhelOS) and another atom y with the thermodynamic 
coupling parameter 5 for the open-conductive (oo) and 
inactivated states (oi) , respectively. The WCA-repulsive 
LJ potential is a very short-range function. For this 
reason, the derivative du^''{r;s) / ds is nonzero only at 
distances that are shorter than or comparable to van 
der Waals contact. It follows that the integrand of Eq. 4 
can contribute to AAGrep only if plr,°'(r;i) is larger than 
pI^'V; s) zero at very short distances. In other words, more 
atoms in the structure must overlap sterically with the 
aromatic side chain of PhelOS in the open-conductive 
(oo) state than in the inactivated (oi) state. This in- 
formation is accessible only through this type of analy- 
sis; it cannot be deduced from the properties of the 



endpoint states (5=0 and 5=1). Searching the con- 
figurations in the FEP/MD trajectories generated 
with 5 = 0.5 reveals that the dominant contributions 
arise from the interaction of PhelOS with Thr75 within 
the same subunit and He 100 in the neighboring sub- 
unit (Fig. 1, B and C). Less important steric contacts 
include Leul05 in the same subunit and Val97 in the 
neighboring subunit. Thus, this analysis of FEP/MD tra- 
jectories supports previous conclusions based on the 
x-ray structure and equilibrium MD simulations (Cuello 
etal., 2010a). 




pH 8.0 

Figure 6. Macroscopic pH-dependent currents (10 s) obtained 
at 100 mV from proteoliposomes containing KcsA-WT (black) 
and KcsA-IlOOA mutant (color). The red and blue traces are rep- 
resentatives of the range of kinetics seen with the mutant. On 
average, the IIOOA mutant inactivates more slowly than WT, with 
time constant varying from ^^15,000 ms, although some records 
display an inactivation constant similar to that of WT. 
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A direct prediction of the analysis of FEP/MD trajec- 
tories is that a possible way to decrease the impact of 
PhelOS on the rate of inactivation might be to substi- 
tute a smaller residue at the position of IlelOO. Interest- 
ingly, macroscopic pH-dependent currents for the 
KcsA-IlOOA mutant indicate that inactivation is indeed 
partially relieved in this channel. The experimentally 
observed effect of mutating an isoleucine at position 
100 to an alanine on macroscopic pH-dependent currents 
is shown in Fig. 6. The macroscopic current response 
upon a pH jump shows that the IIOOA mutant inacti- 
vates more slowly and to a lesser extent than KcsA-WT. 
The mean inactivation time is ^^2.1 s for the WT chan- 
nel and 8.5 s for the IIOOA mutant (not depicted). 
Although there is some variability (Fig. 6, red and blue 
traces) , the slowdown in inactivation is accompanied by 
an increase in steady-state open probability to ^30-40% 
for the mutant compared with 3-5% for KcsA-WT. The 
IIOOA mutant inactivates more slowly than WT on aver- 
age, although there is some variability in value of the 
time constant. These observations are consistent with 
the Phel03-Ilel00 intersubunit interaction noted in 
the analysis of the FEP/MD calculations. Although IlelOO 
is important, this experimental result confirms that 
Phel03 is the critical residue controlling the allosteric 
coupling between the filter and the intracellular gate. 
Although only the most unfavorable Phel03-Ilel00 inter- 
action can be relieved in the IIOOA mutation, it is likely 
that additional (smaller) unfavorable interactions made 
by Phel03 with other residues remain present. The 
latter are all removed in the F103A mutation, which 
may pardy explain why the impact of the IIOOA mutation 
appears to be smaller than that of the F103A mutation. 

Conclusion 

Activation and inactivation gating are allosterically cou- 
pled. The closed conformation of the intracellular gate 
stabilizes the conductive conformation of the filter, 
whereas an open gate favors the inactivated filter. From 
the viewpoint of microscopic reversibility, however, any 
coupling of two structural elements must act in both di- 
rections, implying the possibility of "cross talk" between 
the intracellular gate and the selectivity filter. Thus, 
stabilizing the selectivity filter in its conductive state also 
promotes the closing of the intracellular gate. In KcsA, 
cross talk between the activation gate and the selectivity 
filter was first documented by electron paramagnetic 
resonance spectroscopy (Perozo et al., 1999), where 
opening of the inner gate leads to subde conforma- 
tional changes in the pore helix and outer vestibule. 
These conformational rearrangements, now known to 
relate to C-type inactivation at the selectivity filter, 
are also coupled to the opening of the activation gate 
(Chakrapani et al., 2007a, b). The present study has 
elucidated some of the molecular determinants re- 
sponsible for coupling the state of the intracellular 



activation gate to the nonconductive inactivated con- 
formation of the selectivity filter. Although only a 
single site-directed mutation was analyzed in detail 
using FEP/MD simulations, the results are encouraging 
and suggest that these types of computations could be 
a powerful approach for shedding light on the complex 
networks of residue-residue interactions underlying 
slow inactivation in channels that have been re- 
vealed by previous experimental studies (Panyi and 
Deutsch, 2006). 

Here, we have been mainly dealing with understanding 
the microscopic interactions involved with the thermo- 
dynamic stability/instability of the filter in the closed 
and open states. The mechanism by which the filter 
actually transitions from conductive to collapsed, which 
presumably would require a specification of the kinetic 
rate and the pathway, is beyond the present study. Iden- 
tifying clearly the factors that affect the endpoints of 
such a path, however, is a necessary first step in elucidat- 
ing the mechanism of inactivation. In particular, the 
present FEP/MD simulations are in relatively good 
agreement with experiments, which increases our confi- 
dence in the use of this type of calculations to study 
channel inactivation. Future work is underway to clarify 
the mechanism of recovery from inactivation using 
computational methods for determining the pathways 
connecting the different functional states of the chan- 
nel (Pan et al., 2008; Gan et al., 2009). 
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